function [PrY_1,PrT1_1,PrT2_1]=SIM_Pr

global c1 c2 phi1 phi2 cy gamma1 gamma2 by b1 b2 cgpa phigpa sigegpa sige1 sige2 sigy sigtheta ThetaDraws CutoffDraws GPAerror


theta=ThetaDraws*sigtheta;
GPA = cgpa+phigpa*theta+GPAerror*sigegpa;

PrT1_1 = normcdf((c1+phi1*theta+b1*GPA));
PrT2_1 = normcdf((c2+phi2*theta+b2*GPA));

T1 = (c1+phi1*theta+b1*GPA+CutoffDraws(:,1)*sige1>0)+0;
T2 = (c2+phi2*theta+b2*GPA+CutoffDraws(:,2)*sige2>0)+0;

bias1 = T1 - normcdf(cy+theta+by*GPA);
bias2 = T2 - normcdf(cy+theta+by*GPA);

PrY_1 = normcdf((cy+theta+bias1*gamma1+bias2*gamma2+by*GPA));
Y = ((cy+theta+bias1*gamma1+bias2*gamma2+by*GPA+CutoffDraws(:,3)*sigy)>0)+0;

T1star = c1+phi1*theta+b1*GPA+CutoffDraws(:,1)*sige1;
T2star = c2+phi2*theta+b2*GPA+CutoffDraws(:,2)*sige2;
Ystar = cy+theta+bias1*gamma1+bias2*gamma2+by*GPA+CutoffDraws(:,3)*sigy;






function a = sim_74(param,o)
global theta e12 eg ey;

gamma1=param(1);  gamma2=param(2);
betay=param(3);   cy=param(4);      sigtheta=param(5);    
c1=param(6);      c1d=param(7);     phi1=param(8);    phi1d=param(9);    beta1=param(10);  beta1d=param(11);
c2=param(12);     c2d=param(13);    phi2=param(14);   phi2d=param(15);   beta2=param(16);  beta2d=param(17);
cg=param(18);     cr=param(19);     cm=param(20);
phig=param(21);   phir=param(22);   phim=param(23);
sigg=param(24);   sigm=param(25);   sigr=param(26);
if o ==1 % other race
 %   gamma1 = gamma1+gamma12; gamma2 = gamma2+gamma22;
    c1 = c1+c1d; c2=c2+c2d;
    phi1 = phi1+phi1d; phi2 = phi2+phi2d;
    beta1=beta1+beta1d; beta2=beta2+beta2d;
end;

g=cg+phig*theta+eg;
Tstar = [c1+phi1*theta+beta1*g,c2+phi2*theta+beta2*g]+e12;
T1star = Tstar(:,1); T2star = Tstar(:,2);
T1=(T1star>=0)+0; T2=(T2star>=0)+0;

b1 = T1-normcdf(cy+theta+betay*g);
b2 = T2-normcdf(cy+theta+betay*g);

Ystar = cy+theta+b1*gamma1+b2*gamma2+betay*g+ey;
Y=(Ystar>=0)+0;
PrY = normcdf(cy+theta+b1*gamma1+b2*gamma2+betay*g);

a=[Y, T1, T2, b1, b2,T1star,T2star,Ystar,PrY];



